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Working with a large basis of covariant derivative-based meson interpolating fields we demonstrate 
the feasibility of reliably extracting multiple excited states using a variational method. The study is 
performed on quenched anisotropic lattices with clover quarks at the charm mass. We demonstrate 
how a knowledge of the continuum limit of a lattice interpolating field can give additional spin- 
assignment information, even at a single lattice spacing, via the overlap factors of interpolating 
field and state. Excited state masses are systematically high with respect to quark potential model 
predictions and, where they exist, experimental states. We conclude that this is most likely a result 
of the quenched approximation. 



I. INTRODUCTION 
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In this paper we discuss an application of the variational method for extraction of spectroscopic information from 
. lattice QCD two-point functions. We focus on the use of a large set of interpolating fields allowing access to all 
' J PC with sufficient redundancy to consider several excited states in each channel. As well as the usual extraction 
tIh ■ of mass information, we are also interested in the overlaps of interpolating fields on to states in the spectrum. This 
information is required if one wishes to extract matrix elements between hadron states from three-point or higher 
\ correlators, an example being the radiative transitions considered in We find in this paper that these overlaps 
■ are also useful in aiding the spin-identification of states provided one has established the continuum behavior of the 
CS| ! interpolating fields. 

^vO i Our application of the variational method utilizes the orthogonality of states to extract excited state information. 

This is particularly powerful when a number of excited states are close to being degenerate - conventional multi- 
Nl . , exponential fitting methods are not good at resolving such states. We are able to collect high statistics by the use of a 
C ' small spatial volume ~ (1.2fm) 3 and by confining our study to a large pseudoscalar mass m ps ~ 3 GeV, corresponding 
, to charmonium. This volume is known to be sufficient to house the ground state for most J PC , as illustrated for the 
"charge" radii extracted in [l[ , but we might expect it to be insufficient for radially excited states and those of high 
t". I spin. We will study the volume dependence of a number of correlators and find, somewhat surprisingly, that there 
^ ■ seems to be no meaningful difference between a 1.2 fm box and a 2.4 fm box. 

A resurgence of interest in charmonium physics is underway, spurred by new results from BaBar, Belle and CLEO. 
This is reflected in increased theoretical effort, in particular there have been recent efforts to apply lattice QCD 
techniques to computation of quantities previously the reserve of potential models and sum-rules, such as radiative 
transitions [l[ and two-photon decays Q. In this paper we extract information about the excited state and high spin 
spectrum, which is both interesting in its own right, and vital for concurrent attempts to extract excited state radiative 
transitions. 

We begin in section II with a discussion of the lattice action and other computational details before in section III 
outlining our choice of a broad set of interpolating fields constructed to be irreducible representations of the lattice 
symmetry group at zero momentum and to have simple, known continuum properties. The variational method used 
to describe the data is explained in section IV with focus on ensuring that the solution accurately describes the input 
data. Spectrum results arc then shown in section V, with some discussion of how continuum spin assignments can be 
aided by considering overlap factors. In three appendices we display properties of our interpolating fields and show 
an application of the variational method to toy data. 
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II. COMPUTATIONAL DETAILS 



The computations were performed in the quenched approximation to QCD, using the Chroma software system 0]. 



We employed 1000 configurations on a 12 6 



lattice generated using an anisotropic Wilson gauge action 



x.i>j x,i 



(i) 



defined in terms of simple plaqucttes P flv (x). Here, f3 is the bare coupling, the time component is labeled by 0, the 
spatial indices i,j run from 1 to 3, and £o is the bare anisotropy. We tune £o to the desired renormalized anisotropy 
£ = a s /at = 3 where a s and at are the spatial and temporal lattice spacings, respectively and find £o = 2.464. The 
temporal lattice spacing obtained from the static quark-antiquark potential is a^ 1 = 6.05(1) GeV. This gauge action 
is expected to have 0(a 2 s , a 2 ) discretization errors. A smaller set of 807 configurations on a 24 3 x 48 lattice were used 
to study finite volume effects. 

Anisotropic lattices as applied to charmonium exploit the fact that while the quark mass scale demands a cut-off 
above ~ 1.5 GeV, the internal three-momentum scale is typically much lower, ~ 500 McV. On our lattice, we can 
have both m c a t and \p\a s reasonably small and a spatial length > 1 fm without requiring very many spatial lattice 
sites. Most of the work presented uses one volume, L s « 1.2 fm while a limited comparison is made with a large 
volume L s ks 2.4 fm; previous charmonium spectrum studies indicate that there are no significant finite volume effects 
for lattices of this size or larger [1,0]. 

The quark propagators were computed using an anisotropic version of the O(a) tadpole-improved Clover action 0, 
ID, Q . The Clover action we used has the 

a s and ap = at, we can express this fermion action in terms of dimensionlcss variables tp 



mass" form with spatial fermion tuning. Defining a lattice spacing where 



a-k 

and F^v = a lli a v F l 



I"' 



as 



a t Q = a t m + Wo7o + j- V] W k jk - \ a ^2 a okhk + j- Y] <JkiF t 
so 



Co 



k<l 



(2) 



The factor v is the bare value of the fermion anisotropy. For the field strength tensor F^, we adopt the standard 
0(a 2 s ,a1) clover leaf definition. Here, the "Wilson" operator has the projector property 



where 



V M /(x) 



1 

2a~,, 



U ti (x)f(x + p,a M ) - Ul(x - fia^)f(x - fta M ) 
U fl (x)f{x + pLa^) + UUx - fan) f{x - fan) - 2f(x) 



(3) 

(4) 
(5) 



The projector property ensures that no doubler states are present in the determined mass spectrum. We have tuned 
(mo, v) so as to yield the desired quark mass and such that the speed of light obtained from the meson dispersion 
relations be one, as discussed later. 

We have used the tree-level conditions on c s and Ct as described in Ref. 0]. In particular, we will choose 



v 



1 

Ct= 2 



a s 



1 



u t uj. 



(G) 



where the ratio a t /a s = l/£ is set to the desired renormalized gauge anisotropy. The tadpole factors u s and u t come 
from the fourth-root of the spatial and temporal plaqucttes, respectively, and take the values 0.814 and 0.980. 

The charm mass is determined by tuning the bare quark mass mo non-perturbatively such that the spin average 
of the lowest 5-wave mesons coincides with its experimental value, i.e. (3mj^ + m^ e )/4 = 3.067 GeV., such that 
mo = 0.0401. We tune v non-perturbatively to satisfy the lattice dispersion relation 



Eipf - EiOf 



ajEjpf - alEjQf 
a 2 s P 2 



= 1. 



(7) 



Keeping all other parameters fixed, wc tunc v to satisfy the above relation to within ~ 1% and find v = 0.867. 
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We have used Dirichlct boundary conditions in the temporal direction placing the source five time-slices from the 
wall. All subsequent time-slice plots are with respect to this source position. 

For comparison purposes, we have also used an anisotropic version of the domain-wall fermion (DWF) action [HIIoj| 
with a domain-wall height a±M = 1.7, a fifth dimensional extent L5 = 16, and a quark mass atm q . The kernel of the 
domain action is the same as in Eq. [2] except without the clover terms, i.e., Ct and c s set to 0. As in the clover case, 
the domain- wall quark mass and fermion anisotropy (m q ,v) are chosen to yield the desired S-wave mass and such 
that the speed of light obtained from the meson dispersion relations is one. 

We have computed smeared-local correlators on a subset of the same quenched lattices using both tadpole Clover 
and domain wall fermion quark actions. The effective masses are shown in figure [T] where we see that apart from the 
domain wall fermion oscillations at small times there is essentially perfect agreement between the two actions. Since 
the domain wall action should be O(a) improved automatically we infer that, at least for the ground state spectrum, 
the tadpole Clover is also 0{a) improved to a good approximation. 
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FIG. 1: Effective masses for tadpole clover (red) and domain wall fermions (shades of blue/purple) on same lattices. Note 
masses are shifted for clarity. Relative mass scale is shown. 



In this study we used quenched lattices. In the spectrum of heavy quark states, where unitarity violation due to the 
absence of closed charm quark loops can be neglected, there are likely to be two main effects of this approximation. 
Firstly there are no virtual D mesons in this theory and hence no thresholds for OZI-allowed decay and secondly 
there is the known issue of the incorrect running of the coupling and an associated ambiguity in the setting of the 
lattice scale. In all cases we compute only the connected contributions to correlators - the effects of disconnected 
contributions in charmonium arc believed to be small [111 



III. INTERPOLATING FIELDS 



It is common in meson spectroscopy calculations to focus on the simplest meson interpolating fields, namely the local 
fermion bilincars, ip(x)Tip(x). By appropriately smearing the operator over space, a large overlap on to the ground 
state of a given quantum number can be obtained. These operators are, however, limited in the J = ±+ , l +± , 1 
with which they have overlap. In addition, with only these operators, while one might try various different smearings 
[l3| . it is not possible to explore higher spin states. For these reasons we opt to extend the operator basis by including 
non-local operators, in particular ones constructed from a number of covariant derivatives acting on the quark fields. 
On a discretised lattice covariant derivatives become finite displacements of quark fields connected by links. 

In the continuum one can construct operators which transform in a particular way under the rotation group, giving 
overlap only on to certain spins. On a discrctized lattice the full rotation group is broken down to a smaller group 
of cubic rotations with a limited number of irreducible representations (irrcps). The operators used in this paper are 
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constructed to be both irreps of the lattice rotation group (at zero momentum) and to have definite, simple forms in 
the continuum limit. In this way we have some information about how these operators will behave in the continuum 
limit - we will see that this potentially aides in the interpretation of the lattice data if we are close enough to restoring 
rotational symmetry. 

Our operators are based upon those in [l4j]. The principal differences are that we extend the basis to include 
all possible zero, one and two derivative operators 1 , generalize the derivatives so that the operators have definite 
charge-conjugation at finite momentum and correct a projection operator to ensure orthogonality within a lattice 
irrep. Finally we develop the continuum overlap formalism expressed in appendix lAl 



A. Derivative-based operators 

In [bi ] single derivative operators of the type fofiDjipi were utilized where 1 and 2 signify the flavor. The right- 
facing arrow indicates that the covariant derivative acts only on ipi. For the case of a single flavor we would like the 
operators to be of definite charge- conjugation - it is easy to see using an integration by parts that the above operators 
do not achieve this at finite momentum: 

J (Pxe^^TjD^ix) = - J &3e& 3 %j){x)T i i D j il)(x) + i Pj J Sxe^^^T^x). (8) 

The first term on the right-hand side has the same charge-conjugation as the left-hand side but the second term, 
lacking the derivative, has the opposite C . This can be fixed by replacing the derivative D by the quark-antiquark 
symmetrised D = D — D. The same integration by parts as above will now not yield a term proportional to the 
three-momentum. 

At the two-derivative level, Ref. [l4| construct two combinations, = e^D^D^ and O; = \cijk\Dj D^. In the 
D-type case we ensure definite charge-conjugation at non-zero momentum by redefining Dj = |eyfe| D j D The 
B-type case does not actually require this extension as one can easily see by expressing the pair of derivatives as the 
sum of a commutator and an anticommutator: 

eijkDjDk = <,,^( n,. I), + = e ljk v D k } = -\e iik F^ k . (9) 

Clearly the B-type operator is so named because it corresponds to the chromomagnetic component of the field-strength 

< — > < — > 

tensor. We complete the two- derivative set by adding Ej = Qijk D j D /. (Qijk is a Clebsch-Gordan coefficient as defined 

in the next section) and the Laplacian V 2 = ^ i D iD j . 

In the continuum, Lorcntz symmetry along with parity and charge-conjugation specify the form of the overlap of an 
operator with a state of given J p . In appendix [X] we tabulate the forms, valid for all three-momenta in Minkowski 
space-time, for the operators we use in this paper. 

We approximate covariant derivatives using finite displacements, including the appropriate links. In the results 
presented below we displaced only by one site, displacing by two sites was found to give correlators that did not differ 
considerably. Our implementation takes the form: 

Vjf(x) = £(Uj(x)f(x+fa.) ~ r; ;.r ./,/.;/(,• ./</.. ;) D J f(x)+0(a 2 s ). (10) 



B. Projection onto lattice irreducible representations 

Our method will be much aided by having a basis of operators that transforms irreducibly under the cubic group 
O. The irreducible representations of O, together with their continuum spin contents, are shown in TableU fl5| . The 
operators are constructed as products of gamma matrices and derivatives as outlined above, and therefore we need to 
project these products to their irreducible components. For operators K and L transforming according to the irreps 
A.K and Al respectively, the product M transforms under an clement R of the cubic group as 

M tj = KiLj -» T^iR^T^i^^K^Lj, , (11) 



1 although not all of these are actually used in this calculation 



A 




J 


Ai 


1 


0,4,6,... 


A 2 


1 


3,6,7,... 


E 


2 


2,4,5,... 


Ti 


3 


1,3,4,... 


T 2 


3 


2,3,4,... 



TABLE I: The table shows the single- valued irreducible representations A of the cubic group O, together with their dimensions 
d,A and continuum spin content J. Additional superscripts are employed to denote charge conjugation C and parity P. 



where T A (i?) is the representation matrix for the element R in the irrep A of dimension d^. The essential tool in 
constructing the irreps is the projection formula 

of ae = — £ r^RMR)^ , (12) 
90 

where {Oj} is a basis of operators that is reducible under the cubic group, U(R) is the representation of the rotation 
R on that basis, and {Of} are a set of operators transforming irreducibly under the cubic group. The indices a and 
(3 refer to the rows of the irreducible representation, and the final step in the procedure is to identify a set of linearly 
independent operators {Of~ : i = 1, . . . , d\} as a basis for the irreducible representation. 
Applying the projection formula 1121 to the operator product of eqn. 111! we obtain 

90 £?o 

Identifying a set of d\ linearly independent operators from the set AfA, we obtain the Clebsch-Gordon coefficients for 
obtaining the irreducible representation A in the product of operators transforming according to the irreps Ak and 

Mf = Q ijk Kf«L*?. (14) 

In Ref. the D and B operators are constructed via Clcbsch-Gordan coefficients of Ti(V) g) 2i(V) to transform 
as T2, T\ irreps respectively The complete decomposition of T\ ®T\ is Ai © T\ © T2 © E with the missing A\ being 
the laplacian and the E being the E operator defined in the previous section. The Clebsch-Gordan coefficients in this 
case are Qy fe where all elements are zero except 

Qm = ^; Qi22 = -^; Q211 = -75; ^222 = -^; Q 233 = ^- (15) 

Examining the zero three- momentum lattice irrep projections of the continuum overlaps in appendix ^ we can see 
that certain operators of, say, T\ type can have overlap with higher spins only at finite a with that overlap vanishing 
as a — > 0. The simplest example is the local fermion bilinear ip^ip which, transforming as T± can have overlap with 
continuum spins 1,3,4... at finite a but can only have overlap with spin-1 in the continuum. Several of our operators 
do retain overlap on to higher spins as a — > 0, e.g. the T2 projection of p x D: (p x D)^ 2 = (-ijki>l j\^kim \ D 1 D m i\). The 
non-zero overlaps in the continuum at rest are (0\{p x B)f 2 12 (0, r)) oc |e,-/fe| (0, r) and (0|(p x B)J" 2 13 (0, r)) oc 
£ijk\tkim\ G^ 1 " 1 {0,r) which are linear in the spin-2 and spin-3 Minkowski polarization vectors (p,r), G A " yp {p,r). 
None of our operators retain overlap on to spin 4 or higher in the continuum. The full set of lattice irrep projections 
is displayed in appendix [Bl 



C. Quark and link field smearing 



Constructing covariant derivatives using finite displacements with link fields can potentially produce rather noisy 
correlators. Any individual link is subject to large UV fluctuations over the ensemble of gauge fields. A suitable 
average over links neighboring a given link can help to reduce this fluctuation and providing one does this in a gauge- 
invariant and rotationally symmetric way, the quantum numbers of the operator will not be changed. This subject is 
discussed in detail in [lq . In figure we show the reduction in correlator noise we were able to obtain by judicious 
choice (p = 0.15, N p = 12) of stout-link smearing parameters. 
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FIG. 2: Effective mass (on 300 configurations) for the correlator (p x ' 
Data are labeled by the stout-link smearing parameters, p,N p . 



)ii — (p x B>)ti with link smearing at source and sink. 



Quark bilinear operators that are local or almost local (e.g. single site displacements) tend to have considerable 
overlap with a tower of states. While this can be helpful in a study like this one attempting to consider excited states, 
it is also useful to be able to emphasize the lower states in the spectrum. This can be achieved by smearing the quark 
fields over space with a gauge-invariant cubic approximation to a rotationally symmetric gaussian with any derivatives 
applied subsequently, 

/ 3 \ N 

- 1£) (i + T^m g + U UA,A ) ^ eCT2v2/4 - ^ 

As well as producing improved plateaus, optimizing the quark smearing parameters docs in some cases reduce noise. 
We show in figure [3] two representative cases. Our final parameter selection was a = 4.2, N = 50 for the local, V and 
B operators and a = 5.0, N = 75 for the B operators. 
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FIG. 3: Effective masses (on 300 configurations) for the correlators (p x V)tj — (p X V)tj and (p x B) Tl -(px B) Tl ■ Stout-link 
smearing is applied with p = 0.15, -/V p = 12 and the data are labeled by the parameters of the quark smearing (a, N) which is 
applied only at the source. 



IV. VARIATIONAL METHOD 



Our spectrum results follow from application of a variational method to matrices of correlators. This method takes 
advantage of the orthogonality of state vectors on a basis of interpolating operators. The basic numerical problem to 
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be solved is of generalized eigenvalue type 

C(t)v a = X a (t)C(t )v a . (17) 

In this expression C{t) is the matrix of correlators at timeslice t, i.e. Cjj(t) = (Oi(t)Oj(Q)). We transform all our 
operators to Minkowski space (but retain imaginary time) and in doing so ensure that this matrix is Hcrmitian (more 
detailed discussion of this point can be found in Appendix [B]) . The generalized eigenvectors are orthonormal on the 
metric C(to), u£<7(io)«0 = S a p. The principal correlators X a (t) are the generalized eigenvalues on a given timeslice - 
they can be shown to behave at large times like [13, EH 

X a (t) = e - m °('-*«)(i + ©(e-IMt*-*"))), (18) 

where m a is the mass of a state labeled by a and 8m is the mass gap to the nearest state to a. 

Performing a Cholcsky decomposition on the Hcrmitian matrix C(to) = LL\ where L is lower diagonal, one can 
cast this as a conventional eigenvalue problem 

L- l C(t)l)-\l)v) a = X a (t)(tfv) a , (19) 

where the eigenvectors are u a = (L'v) a . These eigenvectors form a dim(C) x dim(C) unitary matrix that transforms 
our "trial" operators into the optimum linear combinations to overlap with the dim(C) lightest states accessible to 
them. The eigenvectors should be time-independent to the extent that the correlator is saturated by dim(C) states; 
in practice we solve the eigenvalue problem on each timeslice and obtain eigenvectors that can vary with t. 

Clearly from the principal correlators one can extract information about the mass spectrum. There is also useful in- 
formation in the eigenvectors. The conventional spectral decomposition of a two-point correlator (at zero momentum) 
has the form 

c^t) = {o^gm) = E ^r~ e ~ mat > ( 2 °) 

a a 

where the overlap factor, Zf = (0\Oi\a). Examining the form of the generalized eigenvalue equation with substitution 
of the spectral decomposition one finds that 

Z? = (y- 1 )? % /2^e m °* o/2 . (21) 
The inverse of the eigenvector matrix is trivial to compute owing to the orthonormality property V^C{t$)V = I =>• 

v- 1 = ytc(to). 

In practice we work with a finite space of operators and in this case the parameter to plays an important role. The 
eigenvectors are forced by the solution procedure to be orthogonal on the metric C(to) - this will only be a good 
approximation to the true orthogonality (which in the continuum is defined with an infinite number of states and 
operators) if the correlator at to is dominated by the lightest dim(C) states. As such one should choose t large 
enough that you believe the above statement to be true. We demonstrate this using toy data in Appendix [Cl 

In choosing a to value for a given C(t) there are two factors to take into account - the above discussion suggests 
we should push to out to larger values, where the contributions of higher excited states have decayed exponentially; 
however, as we do so we get into a region where the correlator data are typically noisier. As we can sec from equation 
[T9l the noise on C(to) will enter into the solution of the eigenvalue problem at all timcslices and as such we do not 
want to make to too large. We need a criterion to decide upon an optimum value of to - our choice was to define 
a chi-squared-like quantity gauging how well the generalized eigenvalue solution (with time-independent Z values) 
describes the correlators. At a given to we solve the eigenvalue problem to yields masses (from fits to the principal 
correlators - details of the fitting follow later) and Z's. With these in hand we can reconstruct any correlator matrix 
element using equation 1201 A suitable chi-squared-like quantity can be defined as 

1 t max 

X 2 = TTl T7 ; 1— ; r Y Y {CiM - Cl^-mC^itJMCiAt!) - C^ c -(f)), (22) 

±N(N +l)(t max -t )- ±N(N + 3) ^ itt ,t^ +1 3 J 3 

where N = dim(C) and where C is the data correlation matrix for the correlator computed with jackknifc statistics. 

The optimum value of to is chosen to be that which minimizes the chi-squared-like quantity. In fact, since we solve 
the eigenvalue problem on each timeslice we actually get Z(t); we choose to take the Z values (for a given to) from a 
fixed timeslice tz > t such that the chi-squared-like quantity is minimized at this to- Since the Z(t) are reasonably 
flat the chi-squared variation with tz is fairly mild. Insisting that the Z's are time independent for t > t is a reflection 
of the fact that the only time-dependence in the spectral representation [2"01 is in the exponentials. 
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There are some subtleties that arise because we are numerically solving the generalized eigenvalue problem multiple 
times. One is that the phase of the eigenvector for a given state is not a quantity that can be determined with only 
two-point correlator information. From equation [20J it is clear that one could multiply all Zf for a fixed a by the 
same phase and equation [201 would be unchanged. We solve the eigenvalue problem once for each single elimination 
jackknife sample of our ensemble. In general the eigenproblem solver we used does not always obtain the same phase 
for each jackknife sample, the net effect of which is to make the configuration averages of the Z'b appear to be much 
noisier than they really are. Our fix for this is to apply a phase convention that the largest element of the eigenvector 
u = L'v for a given state should be positive and real. 

Another subtlety can arise in the case that two (or more) states are degenerate within the level of the configuration 
fluctuations. In that case the solver may mis-assign a state label on any given jackknife sample (it orders according 
to the size of the principal correlators). This is not a serious problem for the masses since they are, in this case, 
equal within the error bar anyway, but for the eigenvectors it can be troublesome. On different jackknife samples one 
might see the eigenvector flip between two orthogonal choices with the result that the configuration average appears 
to be noisy. We have examined the data in this study and found that this effect does not occur at any meaningful 
level on our ensemble of 1000 configurations. We do, on occasion, see a flip of the eigenvector between neighboring 
timcsliccs, as we show in an example below. It is for this reason that we do not choose to average over timeslices in 
the determination of Z from Z(t). 

The principal correlators are fitted using either a single exponential or a sum of two exponentials with the constraint 
that X(t =to) = 1. The final timeslice used in the fit is chosen by the requirement that the fractional error be below 
10% and the first is selected by maximizing a fit criterion. This fit criterion also decides between the one or two 
exponential hypotheses. 
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FIG. 4: Chi-squared-like parameter as a function of to for the J 4 1 + channel. 

As a concrete example of our solution scheme, consider the A^[ + channel in which we use the following operator 
basis: (7 5 | sm , us , 7°7 5 |sm,us, h x VA 1 |sm,us, P x B Al \ sra ^ s ) 2 . We solve the generalized eigenvalue problem for all t 
between 1 and 10 - the x 2 -like parameter so determined is shown in figure 2J with a clear minimum being observed 
at t = 7. We show in figure [5] the reconstructed diagonal correlators obtained using the solution at to = 7 (and 
tz = 11). For comparison in figure [6] we show a reconstructed diagonal correlator from the solution with to = 2. This 
variational method is only reliable if the eigenvectors of C{t) are orthogonal on the metric of C(to), which only occurs 
if C{to) is saturated by dim(C) states. In the to = 2 case C(to) is not saturated by the eight states available 3 , and 
we subsequently force the eigenvectors of C(t) to be orthogonal on the "wrong" metric, a truncated metric belonging 
to a larger Hilbert space - this shows up at larger times as a poor description of the data. See also appendix [C] where 
this effect is investigated using toy data. 

Returning to the reconstructed correlators with to = 7, we can see the power of this variational method over 
more conventional multi-exponential fits. Considering for example the final correlator plotted, we see that even at 



2 Where sm and us indicate smeared and unsmeared, respectively 

3 which can also be seen in the shortfall at t < 7 in the to = 7 solution, indicating the need for additional states to fully describe the data 
at all t values 
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FIG. 5: Diagonal correlators and variational solution reconstructions (ground state exponential divided out) in the case to — 7. 
Top to bottom operators are 7 s , 7°7 5 , 61 x Vai, p X Ba 1 - Left column smeared, right column unsmeared. Cumulative state 
contributions are shown by the dashed lines. 
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FIG. 6: Example reconstructed diagonal correlator in the case to — 2. 



timcslicc 7 there arc at least 6 states contributing considerably to the correlator. It is unlikely that a multi-exponential 
fitter would converge to a solution with a sum of 6 exponentials fitting over the range 7-32. There are approximate 
degeneracies in the extracted spectrum and in this case the only distinguishing feature of the states are the iTs or 
cquivalcntly the eigenvectors. Without enforcing orthogonality it is hard to see how one would extract meaningful 
information on these degenerate states. 

In figure [7] wc display the fits to the principal correlators of the lightest 6 states. Shown is the effective mass 
although wc remind the reader that the correlator itself is fitted. 

In figure [8] wc display the extracted Z(t) values for the lightest four states for the smeared tfj^ip operator (points 
with very large noise have been removed). The chi-squarcd-likc parameter is minimized at the point tz — 11 which 
is seen to be consistent with the general trend of the data. Between timeslices 19 and 20 we see a possible example 
of the effect mentioned earlier that the eigenvalues for two near-degenerate states (2nd and 3rd) can flip assignment. 

We note that there is not significant curvature for t < to even though in general there can be, corresponding to 
higher mass states having their time-dependence placed in Z(t) by the variational solver. 

V. SPECTRUM RESULTS 

In this study at zero three-momentum we opt to average correlators over diagonal directions, 



O 



T 2 



c 



c E 



= iE^ 1 

i=i 

3 

= i Vc T2 

i=l 
2 

2 / y 11 

1 = 1 



(23) 
(24) 
(25) 



It is easy to show that for the orthogonaliscd lattice irreps constructed these are the only non-zero entries. In 
general if one does not average in this way one will extract approximately degenerate states corresponding to different 
spin-projections (r in Appendix E]) of a single state. 

In the following we will discuss each PC combination separately. 



A. J++ 

In figure [5] wc display the states extracted. The lowest band of states, at around 3500 MeV, can be identified with 
the near-degenerate Xc0,i,2 states. In a potential model interpretation these states have an internal P-wave and are 
split by a small spin-orbit interaction. In principle, given the continuum spin content of the A\, T\, T2, E lattice irreps, 
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FIG. 7: Effective masses of principal correlators and effective masses of fit solutions (one or two exponential fits) 



these states could all belong to a single 4 ++ state. There are a number of reasons for disfavoring such an assignment. 
In the A\ channel we use only the smeared and unsmeared versions of the ipip operator. In the continuum, the lowest 
dimension operator that has overlap with a spin-4 meson has three covariant derivatives - at finite lattice spacing 
such an operator can mix with ^ip, but must be suppressed by three powers of the lattice spacing. We expect that 
this is sufficient suppression relative to 0(a°) overlaps on to spin-0 states that we can neglect it. In fact none of the 
operators used in the Af + , T 2 f+ , E ++ channels have any overlap with 4 ++ in the continuum limit (see Appendix 
IA|) . but the finite a overlap is not always as suppressed as in the point-like case. 

Support beyond simply mass degeneracy for the common 2 ++ assignment of the lightest state in Ti and E comes 
from the extracted Z-values. Consider for example the operator p x V - this has overlap in the continuum limit, at 
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FIG. 8: Z values for the smeared 7 5 operator. Lowest four levels shown. Phases not relevant 
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FIG. 9: J ++ . Left pane: extracted state masses by zero-momentum lattice irrep. Color coding indicates continuum spin 
assignment (black=0, red=l, green=2, blue=3, grey=undetermined). Right pane: Comparison with experimental masses and 
quark potential model masses. Tabulated masses in MeV. 



zero three-momentum, on to a 2 ++ state as follows: 

(0|V'7 i £ | >|2 ++ (d,r)) = Z e 11 (0,r). (26) 

Now, while we might think of the and E irrcps as being independent on a discrete lattice, their particle content is 
clearly related in the continuum limit, e.g. they share a common Z-value: 

(0|(p x V)^ 2 12++ (0,r)) = |e ijfc |(0|V7^^|2++(d,r)) = Z\e ijk \ e> k (0,r) (27) 
(0|(p x V)y2++(0>)) = ® ljk (0\4> 1 W k 4'\2 ++ (0,r)) = ZQ ijk & k (0,r). (28) 

We might reasonably expect that if our simulation can be considered to be "close" to the continuum limit, the Z- values 
extracted from the T2 and E channels would be related as above, up to hopefully small corrections in powers of a. 
This is in fact what we find to a high accuracy. As a result of the correlator direction averaging described above, the 
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relation of the extracted Z's to the Z defined above is Zt 2 = V%Z; Ze = Z, For the lightest state in T2, E we find 
for the smeared operator that 7=^ = 1.00(1). We take this as evidence that these states are rather close to being 
components of the same 2++ state. 

In addition, for the unsmeared p x V operators wc find = 1.00(3). An equivalent analysis can be applied to 

the Oi xl operator yielding 5j = 0.99(3) for smeared and 5^ = 1.12(5) for unsmeared, all of which appear to 
be compatible with a common state assignment. 

The second band of states (around 4100 MeV) features a state in the A2 irrep, which we associate with a spin-3 
state in the continuum (neglecting spin-6 or higher possibilities). This state should be partnered by nearby states in 
the Ti and T2 irreps for which there are candidates being either the first or second excited states in each channel. We 
can use the fact that in the continuum, the Ti and Ai projections of the operator a\ x D share a common Z-overlap 
onto the spin-3 state: 

(0|( ai xD) A2 |3++(0,r)) =|e^|(0|^7 5 DV|3 ++ (0,r)) = Z\e ijk \ G ijk (0, r) (29) 
(0|(oi x D)^ 2 |3++(0,r)) = e^(0|V7 J 7 5 ©V|3 ++ (0,r)) = Ze ljk \e klm \ €? lm (6». (30) 

With the direction averaging used we find that Za 2 = V&Z and Zt 2 = \f\%- The extracted Z-values for the smeared 
operator give 

= 0-15(15); *f*#^ = 0.97(9). (31) 
3Z T2 (l st ex.) v '' 3Z T2 (2 nd ex.) w y ' 

Hence it would appear that the 2 nd excited state in the Ti channel should be partnered with the ground state in Ai 
as components of a spin-3 meson. The remaining components should reside in the T\ channel. In fact none of the 
operators used in that channel have overlap with 3 ++ in the continuum, but overlap can occur at finite a through 
mixing with continuum operators of equal or higher dimension. For example, wc have included an operator {a\ x D) ^ 
which can mix with the same-dimension operator (a± x E)^ which does have overlap with 3 ++ in the continuum. 
The degree of mixing is presumably related to a power of <?(a) with at most logarithmic divergence with a. At this 
stage we will assume that one of the states around 4100 in T\ is a component of this spin-3 state. 

The remaining states in this band then can either be components of a spin-4 meson or a nearly degenerate 
++ ,l ++ ,2 ++ set. As before an argument against the 4 ++ option is that none of the operators used has overlap 
with spin-4 in the continuum - while this is likely to be a strong constraint in the Ai, ipip case, it is less convincing 
in the other channels where we include operators featuring two derivatives. Here the spin-4 overlap may only be 
suppressed by one power of a. Explicitly we find that the overlap of the unsmeared ■ipijj operator onto the first excited 
state in A\ is of the same magnitude as the overlap onto the ground state. 

We can take the alternate route of looking for support for the ++ , 1 ++ , 2 ++ hypothesis. Wc have already seen that 
there should be a definite relation between the Z's in E and T 2 if the state corresponds to components of a spin-2 
meson. From the fact that we have already assigned the 2 nd excited state in the T 2 channel to a spin-3 meson, we'd 
expect that the 1 st excited states in T2 and E would have Z matching. In fact we find 

z z z z 

: (pxV;sm) = 0.69(28); -^-(pxV; us) = 0.95(27); (01 X D; sm) = 0.58(57); -p^_(ai xD; us) = 0.60(24). 



\/2Ze ' V%Ze ' V%Ze ' V%Ze 

(32) 

While not as convincing as the ground state, we do at least see that these numbers are not inconsistent with the 
spin-2 hypothesis. 

Giving spin assignments to higher states becomes increasingly difficult. It is no longer possible to find complete 
sets of components of a given spin across lattice irreps, since in certain channels wc have insufficient operators to 
extract further states (e.g. only two operators in A\ and hence only two states). Additionally both masses and Z's 
get noisier as we go higher in the spectrum. We can tentatively assign the states at around 4300 MeV in T 2 and E to 
a 2++ on the basis of their Z{p x V; us) where = 0.98(17). 

We summarize our assignments in figure [9] where they are compared to experimental charmonium states and the 
potential model states of [19j . Discussion will follow in the next section. 



B. J" 



Results are displayed in figure ITU1 The isolated lightest state in the T\ channel can be identified as the J/ip. In 
this channel, there is then a large gap to the first excited state which is close to two nearly degenerate states. Wc 
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anticipate the explanation of these levels as being two 1 states and a 3 state. In this case we believe that this 
is the first time that a lattice QCD calculation has observed something like the experimental ■0(3686)/V'(377O) pair. 
The extraction of such nearly degenerate states is made possible by the orthogonality properties of the variational 
method. 
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FIG. 10: J . Left pane: extracted state masses by zero-momentum lattice irrep. Color coding indicates continuum spin 
assignment (black=0, red=l, green=2, blue=3, grey,brown=undetermined) . Right pane: Comparison with experimental masses 
and quark potential model masses. Tabulated masses in MeV. 



The band of states around 3800 MeV has a natural explanation as a nearly degenerate set containing a 3 state, 
a 2 state and two 1 states. The first excited state in the channel is associated with the ground state in E as 
a 2 meson using the Z-relation since ^^ _ ( a i x V);sm) = 1.03(13). 

The ground state in A2 can be associated with the ground state in T2 on the basis of mass degeneracy. The 
^-relation is at the three-sigma edge of validity, = 0.76(8). It is quite possible that this discrepancy is a 

discretisation effect, after all the Z-rclations are only supposed to hold precisely in the continuum limit. To complete 
the 3 we associate one of the levels in T\ - since none of the T\ operators we used have overlap with 3 in the 
continuum we cannot apply any Z analysis here, although we do expect to get overlap through mixing with (p x E)^ 
which has overlap with 3 in the continuum. 

A small overlap onto the unsmeared ipjjip operator might be taken as signal for spin-3 nature owing to the sup- 
pression of such an overlap at C(a 3 ), however there are also good continuum physics reasons why a spin-1 meson 
might have a small value of this overlap. In the non-relativistic limit this overlap measures the wavefunction at 
the origin of a quark-antiquark pair - while this can be considerable for an S'-wave state, it is zero for a 3 Z?i state. 
Relativistic corrections convert this zero to a suppression relative to S'-wave states, see for example the experimental 
/^(3686)(2 3 S 1 ) = 279(8)MeV and ^(3770) ( 3 ^i?) = 99(20)MeV. 

We find that the ground state and first excited state in T\ have comparable decay constants (463(8) MeV and 
416(73) MeV) while the second and third excited states have decay constants consistent with zero (182(211) MeV and 
40(153) MeV). This could be explained if they are dominantly 3 D\ and 3 i?3 states. This could be investigated with 
simulations at smaller a where we would expect one / to remain finite (but small) while the other went to zero 4 . 

It is worth noting that these decay constants 5 arc larger than the experimental values. A similar result was seen 



4 provided we included an operator with continuum overlap on to 3 , otherwise the state would decouple altogether 

5 See [lj for details of scale setting of such quantities on an anisotropic lattice 
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for the ground state in the preliminary Np = 2 + 1 results of [2fJ. We are somewhat surprised that this value is 
so high considering the relative experimental agreement using domain wall fermions on the same lattices in [l|. Wc 
might propose that since we are not improving our vector current operator we cannot expect the full O(a) Clover 
improvement, however in the ma s improvement scheme proposed in [2l| there is no improvement at p — (000). 

We have not been able to make spin assignments for the higher states in the spectrum for the reasons outlined 
in the J ++ section above. The only A{~ state we were able to extract is rather close to our temporal cutoff scale 
(and is compatible with the mass reported in [22|) - such a state could be an exotic , but could equally well be a 
non-exotic 4 . We will not comment upon it further. 

We summarize our assignments in figure [10] where they are compared to experimental charmonium states and the 
potential model states of [lj|. Note that we only show those states for which we are confident of the spin assignment, 
the reader may assign the higher states are their own risk. 

C. j+- 

The J" 1 sector has the interesting property that all states in with J-even are exotic in the sense of being inaccessible 
to a fermion-antifermion bound state. Such states can be constructed from higher Fock states and as such are often 
described as being "multiquarks" (extra quark degrees of freedom) or "hybrids" (extra gluonic degrees of freedom). 
Adding an extra pair of charm quarks would take the state mass up to around 6 GeV, which is at the scale of our cutoff 
and where our quenched non-unitarity might be felt. In the physical spectrum it may be possible for light-quarks to 
play a non-trivial role, in this quenched study we can say nothing about this possibility. If a non-trivial gluonic field 
produces exotic quantum numbered states we have hope of seeing it here. 
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FIG. 11: J H . Left pane: extracted state masses by zero-momentum lattice irrep. Color coding indicates continuum spin 
assignment (black=0, red=l, green=2, blue=3, grey=undetermined) . Right pane: Comparison with experimental masses and 
quark potential model masses. Hatching indicates that the spin assignment is not based upon a Z-analysis. Tabulated masses 
in MeV. 



Our results are displayed in figure QT] The lightest state in T\ can be taken to be the l + ~ 
partners in other lattice irreps. We assign the set of states around 4000 MeV to a l +_ . 



?2, A2 pairing comes from analysis of the overlaps where §^r(&i x ©; sm) = 1.07(7). Possible support for the lighter 

of the Ti states being a 3" 1 component comes from the overlap on the ip^ijjip operator which is consistent with zero 
as might be expected for an 0(a 2 ) suppression. The other state has an overlap which is rather similar to the ground 
state and we propose that it is a l -1 . 



h c since it has no obvious 
3 H pair, support for the 
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The states around 4500 MeV are rendered somewhat non-trivial by our operator-limited knowledge of the A2 
channel. There might be other spin-3 states in this channel, but without more operators we cannot tell. If there 
where to be another state near 4500 MeV in A2 we would have two plausible explanations of the states, cither 3 H , 4 H 
or 3 H ,2^ , 1 H , H . We prefer the second choice and have some support for this from the overlaps. The pairing 
shown in blue satisfies Jfo (e&i x B;sm) = 1.03(10) and we identify it with a 2 H . The higher state in Ti have a 

large overlap onto iplilj'tp ancl we propose it is l^ 1 , while the lower is consistent with zero so we tentatively assign 
it to 3 -1 , partnered by the lower state in T2 and a missing state in Ai- The remaining state in A\ is then proposed 
to be a + ~. 

Our spin-assigned spectrum is compared with experiment and quark models in figure 1111 The exotic -1 , 2" 1 
states are somewhat lighter than the signals reported in [3, H3] • 

D. J"+ 

The J v sector again houses exotics, this time in J-odd channels. We can safely assign the lowest two states 
in A\ to '". The lightest states in T 2 ,E appear to make up a 2 '", and Z comparison appears to confirm this: 
-&-(&i x V;sm) = 1.07(8). Above this spin assignment becomes treacherous. We will work under the assumption 

that the unexplored A 2 channel is in fact empty in our mass range - 3 ^ states are exotic and of high spin and might 
reasonably be expected to be heavy. With this we are left with two possible assignments of the states near 4300 
MeV. The first would be to have two nearly degenerate h states, an exotic 1 h state and a 2 h state. The second 
possibility is to have a single h and a 4 h . Analysis of the Z's in T 2 ,E indicates the spin-2 interpretation is at the 
borderline three sigma level: {p x B; sm) = 0.78(7). The overlap of A\ states onto ^7570 is large for the ground 

state, the first excited state and the third excited state, but is consistent with zero for the second excited state. This 
might be evidence for a spin-4 nature for this state, with the overlap suppressed at 0(a 3 ). 

We are not able with the information we have to decisively state which of the two hypothesis above is correct, and 
as such in figure rT2] we show two possible spectrum interpretations of our data. 

In figure [13] we show the effective mass of our smeared (p x -B)t\ correlator along with the correlator taken from 
[lH computed on lattices with the same a s but with different anisotropy (£ = 2). It is notable that the data are 
consistent although the fits are not. 

E. Discussion 

In our results we see thatthe masses of many excited states and states of higher spin are high with respect to 
the quark potential model and, where available, experimental states. We can identify three possible reasons for the 
systematic difference which are the three approximations in our lattice study: finite volume, non-zero a and the 
quenched approximation 6 .. 

The effects of a finite volume of (1.2 fm) 3 would seem to be a likely culprit. Within potential models the wavefunction 
of increasingly excited and higher spin states gets support at larger and larger distances. Estimates from the potential 
model of Barnes, Godfrey & Swanson [23| indicate that for example the rms radius of the first excited XcO is already 
over 1 fm. In order to consider this possibility we computed a limited set of correlators on a quenched 24 3 x 48 lattice 
with all parameters identical to the simulation on 12 3 x 48. In figure [H] we show the ratio of these correlators in 
several typical cases where it is clear that there is no statistically significant difference, particularly in the t > to 
region. To be clear we are unable to see any relevant finite- volume differences between the (1.2 fm) 3 box and the 
(2.4 fm) 3 box. 

We did not study multiple lattice spacings but we can consider results in the literature for the a-dependence of state 



masses. The studies in |24j and [14| use rather similar actions to ours on quenched lattices and see a-dependences that 
are relatively mild if the action is truly improved so that an extrapolation in a 2 s can be performed. If a non-negligible 
0(a s ) term remains the dependence may be more significant. There typically appears to be a decrease in state masses 
as the continuum is approached which may explain part of our high masses. 

In comparing to the quark-potential model the effect of thresholds can be neglected since their effect does not 
appear in this model. Hence we might expect that the principal effect of the quenched approximation would be in 
the incorrect running of the coupling and the associated scale-setting ambiguity. The ambiguity in scale-setting can 



and, in principle, the effect of disconnected contributions to the correlators 
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FIG. 12: J~ + . Top pane: extracted state masses by zero-momentum lattice irrep. Color coding indicates continuum spin 
assignment (black=0, red=l, green=2, blue=3, grey=undetermined) . Lower panes: Comparison with experimental masses and 
quark potential model masses. Hatching indicates that the spin assignment is not based upon a Z-analysis. Lower left pane 
shows the 4 _+ lighter hypothesis. Lower right pane shows the 1 lighter hypothesis. Tabulated masses in MeV 



be clearly seen in the three graphs of figure 13 in [24(. Since the wavefunctions of higher and higher excited states 
sample larger distances we might expect to feel even larger effects in our study. We suspect that quenching may be a 
major contributor to our systematically high masses. 



VI. CONCLUSIONS 



We have investigated the use of a large basis of interpolating fields on the extraction of excited charmonium meson 
states in many J PC channels. We propose an application of a variational method that systematically selects the 
parameter to to ensure the best possible description of the data at all times greater than to. We use the eigenvalues 
(principal correlators) to determine the mass spectrum and the eigenvectors to determine the overlaps of our inter- 
polating fields on to the states extracted. These numbers arc used herein to aid continuum spin assignment and are 
required quantities for study of three-point functions with the aim of extracting transition matrix elements. 
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FIG. 13: Effective mass of (p x B)t± correlator from this analysis and from [lj- Upper horizontal line is result of fit performed 
in [3. 



The variational method has the important feature of using the orthogonality in a space of interpolating fields of 
different states - this is a powerful lever-arm in the extraction of near-degenerate states. In principle the method 
can miss states that appear in the spectrum if no linear combination of the interpolating fields used has a sufficient 
overlap on to the state in question; however, this is equally true of any "fitting" method and simply suggests the use of 
the broadest set of interpolating fields possible. With perfect statistics this analysis method should not produce any 
"false" , additional states in the spectrum as can happen in simple non-linear exponential fitting schemes where one 
adds states to minimize a chi-squared. With finite statistics and a simple ordering of eigenvalues there is the possibility 
that configuration-by-configuration there is state mis-assignment which can appear as several levels appearing to be 
degenerate within large errors. We did not observe effects like this on our sample of 1000 configurations, but with 
a smaller ensemble and the consequent larger fluctuations of the jackknife bins it might be expected to occur. The 
effect would be partnered by large fluctuations, configuration-by-configuration, of the eigenvectors, and a possible 
mechanism to control it would be to order the levels not simply by the size of the eigenvalue but by the similarity 
(via a dot-product, say) to some established basis of eigenvectors (for example comparing each single-elimination bin 
to the eigenvectors on the entire ensemble average). 

Our interpolating fields are designed so that they have relatively simple and known overlap structures in the 
continuum limit, this was used to aid in the continuum spin identification of states. Supplemented with analysis at 
more than one lattice spacing we believe this will remove much of the possible spin ambiguity of working on a cubic 
lattice. 

We observe in our results that our excited state masses are systematically high with respect to experiment (where 
it is measured) and to the potential models of fl9l ]. Our best guess is that this is a combination of not having 
extrapolated to the continuum limit and the effect of the quenched approximation. Both these issues can be remedied 
with further computation now that the efficacy of the fitting model has been demonstrated. 

Arguably the most intriguing channel considered in this study i s l ". The lowest spin contributing to this is the 
exotic 1 h which has been the subject of several analyses [22l. [2o1. [26T l27l. |28| . In this study we found that we could 
not determine whether the lightest state in this channel (at ~ 4300 MeV) is exotic spin-1 or non-exotic spin-4. In the 
potential models of [l{|, a non-exotic 4 ^ state is expected in this mass region. With the uncontrolled systematics 
related to quenching in this work we can have only limited bearing on phenomenology, but we lay out here the 
potential for misinterpretation in this channel and look forward to a relevant discussion of these matters in a future 
analysis with improved systematics. We note that an analysis varying the lattice spacing using the same operator set 
might shed light on the spin-assignment - no operators are used with have overlap on to spin-4 in the continuum, any 
overlap with such a state is suppressed by powers of a that will reduce as a — > 0. Thus a rapidly reducing Z(a) would 
suggest a spin-4 interpretation. 
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(p x V)t 2 , smeared-unsmeared (ai x V)t 2 - Also shown are the optimum to found in the 12 3 analysis. 
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APPENDIX A: CONTINUUM OVERLAPS 



In this appendix we tabulate the Lorentz covariant kinematic structures corresponding to overlap of fcrmion bilinear 
interpolating fields with a limited number of Lorentz indices and states of definite J PC . In all cases the fermion fields 
should be considered to be evaluated at the origin of Minkowski space time. 



1. No Lorentz Indices 



(0\^\0++(p)) = z 
(0\^ 5 4>\0~ + (p)) = z 
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2. One Lorentz index 



(o|v;o^v|o pc (p)) = zp? 

(0|^V|i (_p)c (p,7-)> = Z £ p (p,r) 



3. Two Lorentz indices 



(0|^^V|0 PC (p)> = Z g^ + Z vV » V v 
(O|^0^|l rc (p>)) = Ze^ afJ e a (p,r)p p 
(0|^"^|2 PC (p>)) =Ze^ (p,r) 
(0\i>O^iP\l^ c (p,r)) = Z + { G"(p>K+ e v {p,r)p») + Z_( ^{p,r)p v - 

[y^y] 7 ^£)^ 7 5 T tl Z) 1 ' {£>\ Z) 1 '} 7 5 {L>^, I? 1 '} j 5 (I)p, Z?^] 7 5 [L^\ D"] 
+- ++ ++ -+ +- 



4. Three Lorentz indices 



(0|^O^^|0 PC (p)) = Z Q g^pP + Z+{g^p v + ^V) + Z.{g^p v - g vp pP) + Z p pW 
(0|i/iC" Jl/ ^|0 ( ~ P)C (p)) = Ze» vpa p a 

(0\$O^\l PC (p,r)) = Z e^" 6 (?,r) + 6^»/ 

+ Z+(e^f> a e/3(p,rK + e" pQ V ^( P >K) + Z-(^p a &p{p,r)p v - e v P^p a e f3 {p,r)pP) 
(0|^O^^V|1 ( " P)C (P, r)) = V G p (p, r) + G"(p, r) + g vp e"(£ r)) + Z_ (.g^ e'(S r) - G^(p, r)) 

+ ZpapV e p (p, r) + Z p+ (pV e"(fl r) + pV e M (£ r)) + ( P y e"(£ r) - P >" g p (p, r)) 
(O|^0^|2 PC (p>)) = Z Q e^(p,r)p" + Z + ( & ^(f, r )p u + e up (p,r)p») + Z_( £"'{?»/- e" p (p»/) 
(0|^"V|2 ( - P)C (P, r)) = Z e^^ e>(jf, r) P/9 + Z+(e^ r) P/J + e v "^ G£(p, r) P/9 ) 

+ Z_ (e^ G* (p, r) P/3 - e"^ G£(p, r)P0) 
(0|^C^|3 ( - p)c (p,r)) = Z e^ p (p,r) 



Q^vp 



DP {DP, D V }Y {DP, D"}j 5 Y [DP, D v YfP [DP, D v ]j 5 Y 



++ +- -+ ++ 

5. Four Lorentz indices 

Here we limit ourselves to the case where the first two indices are antisymmetric and the final two have definite 
symmetry. 
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(0\^ p ^\Q rc (p)) = Z (W - 9 VP D + Z P± [{g pp p"p a - g v W) ± (<f<W - <TW) 
(0|^C" il/ '" T i/;|0 ( - p)c (p)) = Zoe^ 1 ^ + Z p {e^ paa p a p v - e up<ra p a p*) + Z±{f vpa p a p a ± e^" ro, p ot p p ) 
(0|^O^V|l PC (p; r)) - Z p± (e^« e a (fl r)p CT ± e^ CTQ 6„(?, r)p>) + Z e± (e^ pa p a E°{p, r) ± e^ CT > Q e p (p, r)) 
+ e a (p, r)p v ± 6 a ® r)/) + Z e (e^> Q r) ± e"™> Q e"(p, r)) 

+ Z g g p °^ £ a (p, r)p p + Z p *p p p°e^ £ Q (p, r)p p 
+ Z gp± [(^ e a (p,r) P/3 - .g^e^ G Q (p,r) P/3 ) ± (g^ e » pa P e a (p,r)p p - g^e^ e a (p,r) P p) 

+ Z P 2 ± [(//e^ e„( P »|, rP Y f ^ e Q (p,r)p^) ± {fp°z vp <*V 6a(p»W-pY^ r) P/3 ) 
(0|^^V|l ( - P)C (p,r-)) = V( e"{p»p"- + Z p2 p p p°( e p (p,r)p»- £ v (p,r)p») 

+ Z 1± [{g pp E»(p, r)p° - g vp E p (p, r)p°) ± {gT € v & r) P p - g va r)p") 
+ ^± G P (P, r)p" - g p ° €"(?, r)p v ) ± (g vp e CT (£ r)p p - <f» r)/) 

(0|ViO^V|2 PC (p, r)) = Z 9± r) - g vp r)) ± (5^ e""^, r) - ^ e w g r)) 

+ Z p2± [(pV G w (p; r) - P y e"*® r )) ± (p p p a e" p (p, r) - pV r)) 
(0|^" CT V|2 ( - P)C (p,r)) = Z ± (e pvpa E° a {p,r) ± e p (p») + Z a (e ppaa G») - e£(p») 



Zi± 



(0\^O p » p ^\3 PV (p,r)) = Z e p » af3 p a ef (?,r) 



Z± 



{e ppa ^ Pa € v f{p,r) - e vp ^p a eg*(jf,r)) ± (e p ™V e*"(p» - e v ^p a €% p (?,r)) 

(0|VSO^V|3 ( - p)C (p;r)) = Zod""*?^ e vpa (p,r)-p' y E ppa (p,r)) 
(0|^O^ p,T V|4 rC (p,r)) = 



Q^vpa 



j p Y[D p , D a \ l p Y{D p 1 D a } 
+- ++ 



6. Minkowski hermitian operators 

The following operators are hermitian in Minkowski space. Using these operators one is guaranteed a hermitian 
correlator matrix which is a requirement for the variational method. 

ipip ipij 5 ip Tp^f p ip ipj 5 j p ip ipij^j^ip. 

tpUiiJ) Tpij^ip ^Dtv ^7 5 7 M Bl'0 ^7 m 7 1/ SIV'- 
ipElip i/iz7 5 Elv rpy^ip V'7 5 7 M El'0 ^i^^Elip. 



APPENDIX B: LATTICE IRREP OPERATORS AT p— (000) 
Third column indicates the quantum numbers accessible in the continuum limit. 



op. 


name 


cont. 


op. 


V* 


(a x V)ti 


1— 


D* 


7 5 V* 


(tt x V) Tl 


1+- 


7 5 D l 


7 4 7 5 V i 


(tT( 2 ) X V)ti 


1+- 


7 4 7 5 D l 


7 4 V ; 


(a ( 2 ) x V)t, 


1-+ 


7 4 D i 



name 



cont. 


op. 


name 


cont. 


2++ 


M i 


(a x B) Tl 


1+- 


2"+ 


7 5 B' 


(tt x B) Tl 


1— 


2-+ 


7 4 7 5 B i 


(7T( 2 ) X B) Tl 


1~ 


2+- 


7 4 B' 


(a 0(2 ) x B) Tl 


1++ 



7*V l 
e 4 ife7 J V fc 
|e iifc |7 J 'V fc 

7 4 7 i V i 
e^7V'V fc 
|e yfc |7V'V fe 
Q^7V'V fc 

7 5 7 i V i 
e tJ fc7 5 7- , V A; 
|e yfc |7 5 7 J V fc 
Q«fe7 5 7 j V & 



(p x V) Al 
(P x V) Tl 
(p x V) T2 
(pxV) £ 

(P(2) X W) Al 
(P(2) X V) Tl 
(P(2) X V)t 2 
(P(2) X W) E 



0++ 




1++ 




2++ 


e yfc 7 3 'E) fc 


2++ 




0++ 


7 4 7 J D* 


1++ 





2++ 
2++ 



(ai x V) Al 0' 

(ai x V) Tl 1" 

(ai x V)t 2 2" 

(ai x V) B 2" 



(a x D)t 2 
(tt x D) T2 

(7T (2 ) X D)t 2 

(a 0(2 ) x D) T2 



(p x P) A2 3~ 

(p x D) Tl I" 

(pxB) T2 (2,3)" 

(pxD)js 2— 

(p (2 ) x D)a 2 3"" 

(p (2 ) x D) Tl 1"" 

(p (2) xD) T2 (2,3)" 

(p(2) x B) E 2— 



|e ljfc | 7 5 7-''lD) fe 
Qy fc 7 5 7 J D fc 



7 4 7 5 7*V* 


(h x W) Al 


o-+ 


e l ^7 4 7 5 7 J V fc 


(h x V) Tl 


1-+ 


|eiifc|7 4 7 5 7 j V ft 


(h x V) T2 


2-+ 


Q yfc 7 4 7 5 7 J V fc 


(&i x V) B 


2-+ 


op. 


name 


cont. 


E* 


(a x E) Ta 


2++ 


7 5 E l 


(tt x E) T2 


2-+ 


7 4 7 5 E 4 


(tT( 2 ) x E) T2 


2-+ 


7 4 E l 


(a 0(2 ) x E) T2 


2+- 




(p x E) Tl 


(1.3)- 




(p x E) T2 


(2,3)" 


R y - fc7 VE fe 


(p (2 ) x E) Tl 


(1.3)" 


T ijk ^E k 


(p (2 ) X E) T2 


(2,3)- 



7 4 7 5 7 i D i 
|£yfe|7 4 7 5 7 J D 



R tjk j 5 1 m k (ai x E) Tl 
T ljfc 7 5 7^E fc (a x x E) T2 



^■ fc 7 4 7 V'E fc (ME) 



1,3)++ 
2,3)++ 

1,3)+" 
2,3)+" 



(oi x B) A2 
(ai x B) Tl 
(oi x D) T2 
(oi x D) E 



(h xl 

fax: 

(&1 x 



0a 2 

»)t 2 
5)b 



3++ 
1++ 
(2,3)++ 
2++ 

3+- 
1+- 
(2,3)+" 
2+- 



7"B 
eiife7 j 



7 4 7 l 



ey fc 7 4 7 J 'B fc 
|e Jjfc |7VB fe 
Qy fc 7 4 " JMA; 



7 5 7 i B i 

eiife7V'B fe 
|e JJ fc|7 5 7 J B fe 



7 4 7 5 7 i B i 
eyfe7 4 7 5 7 J B fc 



Jyfc7 4 7 5 7 J I 



(pxB)^ 0"+ 

(p x B) Tl 1"+ 

(p x B) T2 2-+ 

(p x B) B 2-+ 

(p (2) x B) Aa 0"+ 

(p ( 2)XB) Tl 1-+ 

(p (2) x B) T2 2 

(p (2 ) X M) E 



{h x M) Al 
(6i x B) Tl 
(6i x B) T2 
(6i x B) B 



T„fc 7 4 7 5 7 J IE fc (6i x E) T2 
KjjfcjTijfc are the Clebsch-Gordan coefficients for T\ ® E = Ti © T2- We find the non-zero values 



+ 
2-+ 



(oi x B) Al 0+- 

(oi x B) Tl 1+- 

(oi x B) Ta 2+- 

(aixB) E 2+" 



0++ 
1++ 
2++ 
2++ 



> _ VS. 
mi — —, 



»412 



» _ \/3. 
^221 — - — , 



^222 



1. 
4 ' 



^332 



-Hill — 4, Jl 112 — —7 J1 221 — 4, 222 — 4- 1 -1^332 — 
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APPENDIX C: VARIATIONAL ANALYSIS OF TOY DATA 



In this section we consider the variational method applied to a finite spectrum of ten states when using a set of 
only five operators. 

We use the mass spectrum m a =(0.5, 0.6, 0.63, 0.7, 0.72, 0.8, 0.83, 0.9, 0.91, 1.04) proposing that there exist 
ideal operators which have unit-normalized overlap on to only one state, i.e. there are Zf = <5 ?ct . We build a model 
analogous to a set of 'trial' operators (linear combinations of the ideal operators) by multiplying Z by a random 
orthogonal 10 x 10 matrix, Zf = ZJMf. We construct a non-diagonal correlator matrix with elements 



Cij(t) 



2m a 



(CI) 



Application of the variational method to the full 10 x 10 matrix of correlators solves the problem exactly for any 
value of to. The situation we deal with in practice is not like this - we do not have anything like a 'complete' basis of 
operators, we can model this by truncating the correlator matrix down to a 5 x 5 submatrix. 

With a selection of random orthogonal rotations (M|) we observe a range of behaviors in the solution to the 
variational problem, here we show two such cases. 

1. Case 1 

In figure [TBI we show the effective masses of the principal correlators for the choices to = 2, 10, 30. We see that there 
is relatively little sensitivity in this case to the value of to. On the other hand we see considerable sensitivity to to in 
the Z(t) values shown in figure [TBI and only approach accurate determination of the Z for the largest to value. 
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FIG. 15: Extracted meff A(t) for to — 2(red), to = 10(green) and to — 30(blue). Horizontal black lines are the input spectrum. 



Case 2 



In figurefTTlwe show the effective masses of the principal correlators for the choices to = 2, 10, 30, note the qualitative 
differences with respect to the previous case including a flip in level ordering and the lack of a plateau in the 4th 
excited state. The plateaus of the first and second excited states clearly improve as to is increased. 
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FIG. 17: Extracted meff X(t) for to = 2(red), to = lO(green) and to = 30(blue). Horizontal black lines are the input spectrum. 
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